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Coarsening of sand ripples is studied in a one-dimensional stochastic model, where neighboring 
ripples exchange mass with algebraic rates, T(m) ~ m 1 , and ripples of zero mass are removed from 
the system. For 7 < ripples vanish through rare fluctuations and the average ripples mass grows 
as (m)(t) ~ — 7 _1 ln(t). Temporal correlations decay as t~ x ^ 2 or t~ 2 ^ 3 depending on the symmetry 
of the mass transfer, and asymptotically the system is characterized by a product measure. The 
stationary ripple mass distribution is obtained exactly. For 7 > ripple evolution is linearly unstable, 
and the noise in the dynamics is irrelevant. For 7 = 1 the problem is solved on the mean field level, 
but the mean-field theory does not adequately describe the full behavior of the coarsening. In 
particular, it fails to account for the numerically observed universality with respect to the initial 
ripple size distribution. The results are not restricted to sand ripple evolution since the model can 
be mapped to zero range processes, urn models, exclusion processes, and cluster-cluster aggregation. 

PACS numbers: 45.70.Qj, 47.54.+r, 05.45.-a, 05.10.-a 



I. INTRODUCTION AND MOTIVATION 

When a surface of sand is exposed to wind or water 
flow, patterns like ripples or dunes are commonly formed. 
The physics of this process is extremely complex because 
it involves the interaction of a granular medium with a 
possibly turbulent hydrodynamic flow fl. It is therefore 
desirable to develop simplified models that capture some 
of the key features of the pattern formation. 

In this paper we are concerned with a class of mod- 
els which focus on the role of the mass transfer in the 
evolution of the pattern. Along a one-dimensional cut 
perpendicular to the ripples, the pattern is described by 
a set {Xi} of ripple lengths, where the index i labels the 
ripples in the array. The Xi are used here as a general 
measure of ripple size, without reference to the detailed 
geometry of individual ripples (see Fig.|l|). In particular, 
we do not distinguish between the linear size of a ripple 
and the mass it contains (for further discussion of this 
point see ||). 

During the evolution of the patterns, the flow transfers 
mass between neighboring ripples. The central assump- 



FIG. 1: Experimental image of vortex ripples in a one- 
dimensional annular geometry |y . The amplitude of the fluid 
oscillations is denoted by a. The line above the pattern shows 
a fit of triangles with a constant slope. Courtesy of K.H. An- 
dersen. 



tion of the model is that the mass transferred to ripple 
i from ripple i + 1 or i — 1 (per unit time) is a function 
r(Aj) of the size of the ripple which gains the mass. Fur- 
ther motivation for this assumption will be given below. 
We refer to T(X) as the robber function pq ]. 

Depending on the characteristics of the flow, the mass 
transfer between ripples can be symmetric or asymmet- 
ric. In the symmetric case the balance between loss and 
gain processes for a given ripple leads to the evolution 
equation (H 

^ i = 7 i [-r(A l - 1 ) + 2r(A l )-r(A !+1 )], (1) 

while in the asymmetric case (assuming, say, mass trans- 
fer only to the left) one has 



dAi 
dt 



-r(Ai_i) + r(Ai). 



(2) 
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The factor 1/2 in Eq. ([!]) makes the time scales for both 
the dynamics equal. 

A homogeneous state of equally sized ripples, Xi = A, 
is stationary under (Q) and (|^) , but its stability depends 
on the derivative of the robber function: The pattern is 
stable for T'(X) < and unstable for F(A) > |. In 
the unstable case the dominant mode is a modulation of 
period 2, in which every second ripple grows and every 
second one shrinks. As the size of the shrinking ripples 
reaches zero in a finite time, the evolution equations (|l|) 
and (||) have to be supplemented by an extinction rule: 
When the size of a ripple vanishes, it is removed from 
the system and the remaining ripples are relabeled such 
that the previous neighbors of the removed ripple become 
neighbors of each other. Extinction events contribute 
to the coarsening of the pattern, i.e., to an increase of 
the mean wavelength. In this work the reverse process 
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of ripple creation is not considered, hence coarsening is 
irreversible. 

The symmetric mass transfer model (|l|) was first pro- 
posed as a description of vortex ripples in coastal waters, 
which are created under the oscillatory flow of surface 
waves ||. In that context the dependence of the rob- 
ber function on the size of the gaining ripple is moti- 
vated by the observation that the mass transfer is effec- 
tuated mostly by a separation vortex which appears in 
the wake of that ripple. Numerical simulations || and 
experiments [|J show that T(A) is nonmonotonic, with a 
maximum near A = a, where a is the amplitude of the 
fluid oscillations. Thus patterns of wavelength A < a 
(A > a) are unstable (stable), and the main interest is 
in the wavelength selection process starting from a short 
wavelength, unstable state g, Q . 

A related, asymmetric mass transfer model for wind 
driven sand ripples was introduced in The basic hy- 
pothesis of the model is that wind ripples wander with 
a speed that is inversely proportional to their size. This 
implies that a leading ripple (ripple i + 1) is eroded by 
the trailing ripple (ripple i) at a rate which is propor- 
tional to I/A;, so the resulting evolution equation is of 
the type @ with T(A) - I/A. Since T'(A) < 0, the homo- 
geneous pattern is stable. However, when fluctuations are 
included by discretizing the ripple sizes and implement- 
ing a stochastic mass transfer rule, a fluctuation-driven 
coarsening mechanism becomes effective and leads to an 
increase of the mean wavelength with time t as mi. 

In this article we consider a class of stochastic models 
whose noiseless counterparts are described by (Q) or (|^). 
We concentrate on monotonic, algebraic robber functions 
r(A) ~ A 7 and study the coarsening process regarding 7 
as a variable parameter. For 7 < this extends the re- 
sults of H on fluctuation-driven coarsening. The case 
7 > is a simple realization of linearly unstable ripple 
evolution, and it is studied here as a first step towards a 
better understanding of models with nonmonotonic rob- 
ber functions [|[ ||, |j] . Although the models are defined 
using the terminology of sand ripples, they are connected 
to other problems in nonequilibrium statistical physics. 
For example, for 7 = the system maps to coalescing 
random walks and is therefore exactly solvable. Other 
equivalences include exclusion processes, zero range pro- 
cesses, urn models, and cluster-cluster aggregation. 

Our main results are the following. In general, one can 
identify two time scales in the dynamics: The one of rip- 
ple extinctions and the other at which the system would 
equilibrate to a steady state in the absence of extinctions. 
For 7 < the loss of a ripple is a rare fluctuation when the 
mean ripple size is large. Therefore the two time scales 
are well separated, and the system has time to relax to a 
quasi-steady state between ripple extinctions. We show 
that this state is characterized by a product measure. 
This justifies the mean field assumption made in 0], and 
allows us to calculate the stationary ripple size distribu- 
tion. The product measure becomes exact only at the 
limit t — > 00 as the correlations in the system decay as a 



power law. The average ripple size grows logarithmically 
at late times, with a prefactor — 7 . 

For 7 > extinctions are frequent events which occur 
on the same time scale as the evolution of the surviving 
ripples. We find in this case that the noise is irrelevant, 
so that the dynamics can be described by (|l|) and (j^). 
For < 7 < 1 the mean ripple size grows algebraically 
with the exponent 1/(1 — 7), while the growth is exponen- 
tial for 7 = 1. In the latter case the evolution equations 
become linear, and the problem can be solved exactly 
on the mean field level. The mean field theory repro- 
duces the exponential growth for the mean ripple size, 
but incorrectly predicts a dependence of the ripple size 
distribution and the coarsening law on the initial condi- 
tions. 

In the next section the model is introduced and its re- 
lations to other models are discussed. Algebraically de- 



caying robber functions (7 < 0) are considered in Sec. Ill 



The produ ct form of the mass distribution is deri ved in 
Sec. [II A , the coarsening law is calculated in Sec. Ill B , 
and the a pproach to the product measure is analyzed in 
Sec. Ill C| . Section |y| is devoted to algebraically growing 
robber functions (7 > 0). The mean- field theory is first 
developed for 7 = 1 and then compared to simulations 
(Sec. |fVA|). Section ^VB| examines the case < 7 < 1. 
Conclusions and open questions are formulated in Sec- 
tion [yj 



II. THE STOCHASTIC RIPPLE MODEL 
A. Definition and simulation algorithm 

In the stochastic model a sand ripple is characterized 
by its mass m. The mass variables are integers such that 
each ripple consists of mt elementary mass units and oc- 
cupies a site i on a one-dimensional lattice. The mass 
is conserved, i.e., M :— X)j TO i = const.. The rrii corre- 
spond to the length variables A; used in equations ([l]) and 
(^) . As mentioned in the Introduction, the mass and the 
length of ripples are here considered to be indistinguish- 
able. We use different symbols for two reasons. We want 
to make a clear distinction between (i) the real and inte- 
ger valued ripple sizes and (ii) between the deterministic 
and noisy dynamics. 

Ripples interact only by exchanging mass with their 
nearest neighbors with an algebraic mass transfer rates 
r(m) = Tom 7 . Since the constant To affects only the 
time scale it will be set equal to unity from now on. If 
ripples obtain mass only from one of their neighbors, say, 
from the right one, the mass transfer is called (totally) 
asymmetric. If the mass comes from both neighbors we 
call the dynamics symmetric. As was discussed in Sec- 
tion I, the asymmetric mass transfer naturally arises in 
the case of wind ripple formation H whereas the sym- 
metric dynamics takes place for ripple patterns forming 
under an oscillatory flow [|| . 

In addition the model includes the removal of ripples 
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when their mass becomes zero. This is done such that 
lattice sites containing no mass are eliminated from the 
system. In this way each ripple always has a neighbor 
from which it can gain mass. If we denote the number of 
lattice sites at time t by N(t), the average ripple mass is 
(m)(t) = M/N(t). 

In the simulations three different initial conditions are 
used. As random initial conditions we denote the case in 
which the probability to have a ripple of size m is given by 
the geometric distribution (1 — q)q rn ~ 1 , with < q < 1. 
The probability q is related to the mean ripple size as 
(m) = (1 — q) . A distribution m, = (m) Vi is referred 
to as monodisperse. The third possibility is a Poisson 
distribution. 

The dynamics is implemented as follows. First a rip- 
ple is selected randomly and time is incremented by 
N(t)~ r~^ x! where r max is the maximum of all the rates 
of the ripples in the system at time t. Denote the mass 
of the selected ripple by m. If x < r(m)/r max , where 
x is an uniformly distributed random number in the in- 
terval [0, 1], the ripple gets a unit mass from its nearest 
neighbor. Otherwise a new ripple is selected and the pro- 
cess is repeated. For symmetric dynamics the neighbor 
is selected randomly whereas in the asymmetric case it 
is always the right one. 



B. Relation to other models 

The model defined above is inspired by the worm 
model originally introduced to describe the coarsening 
of wind ripples p|. Here it is generalized in two respects. 
First, the mass transfer rate in the worm model is in- 
versely proportional to the ripple mass whereas in the 
generalized model it is given by T(m) ~ to 7 . Second, 
we consider also the symmetric mass transfer between 
neighbors. In the case of asymmetric dynamics and for 
-) = — 1 our model reduces to the worm model. 

Apart from the extinction step, the sand ripple model 
is similar to a zero range process §,01- Both models 
are defined in terms of conserved, integer mass variables 
rrii which interact through the (symmetric or asymmet- 
ric) exchange of unit masses between nearest neighbor 
sites of a lattice. The key difference is that in a zero range 
process the mass transfer rate is a function of the mass 
at the site of departure, while in the sand ripple model it 
depends on the mass at the target site. This reverses the 
sign of the right hand sides of (jf]) and and hence the 
stability properties of the model: In a zero range process 
the homogeneous state is stable if r'(m) > and unsta- 
ble if V {fa) < 0. The coarsening behavior in zero range 
processes with nonmonotonic robber function is relevant 
to clustering in granular gases || . 

The occurrence of irreversible extinction events in our 
model is reminiscent of certain urn models that have 
been proposed in the context of glassy dynamics [ |l0| . 
For example, consider the backgammon model 11 which 
is defined by M particles distributed among N boxes 
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FIG. 2: Mapping between the asymmetric worm model and 
the exclusion process with coalescence. 



with mi particles in the ith box. The Hamiltonian is 
H = — Y^i^mtfi so that the energy corresponds to mi- 
nus the number of empty boxes. The TV-fold degen- 
erate ground state therefore consists of a condensate, 
where all the particles belong to one box. Associating 
the masses of ripples with the particle numbers and lat- 
tice sites with boxes, the ripple evolution becomes similar 
to the backgammon model at zero temperature, where 
the empty boxes are not refilled once they have become 
empty (ripple extinction). 

In contrast to the ripple model, the urn models have no 
spatial structure, i.e., mass transfer is possible between 
any pair of boxes. In the standard dynamical scheme, 
originally due to Ehrenfest, in each time step one of the 
balls is chosen at random and a move to another box is 
attempted fl(| . The probability for a box to be chosen 
is then proportional to its occupation number. In our 
setting this corresponds to a mass transfer rate T(mi) ~ 
rrii, where i is the site of departure; in this respect the urn 
models are related to zero range processes. Since V > 0, 
the homogeneous state is linearly stable and coarsening 
(i.e., evolution towards the ground state) is very slow. 

The sand ripple model can also be mapped to an ex- 
clusion process J?], |l2[ [l3| . The mapping can be done 
in two ways which differ in how the disappearance of 
ripples is taken into account. The mappings proceed 
along the lines of [fl4| and the first one is schematically 
presented in Fig. H One constructs a new lattice with 
L(t) = M + N(t) sites. The mass variables of the rip- 
ples turn to rrii consecutive holes separated by particles 
on the new lattice. More precisely, there exist particles 
on sites i + 531=1 m k (* = lj • • • ; AT(i)) while the rest 
are empty. Moving one mass unit from one ripple to its 
neighbor corresponds to a hop of a particle in the ex- 
clusion process. Naturally the exclusion process is either 
symmetric or asymmetric as the original dynamics. For 
the exclusion process the loss of a ripple becomes a coa- 
lescence of particles at contact, which changes the length 
L(t) of the system. 

As the masses map to holes, the hopping rates of parti- 
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FIG. 3: Mapping between the worm model, the asymmetric 
exclusion process (ASEP), and the cluster-cluster aggregation 
(CCA). 



cles depend on the distances between them. For negative 
(positive) 7 there is a repulsive (attractive) interaction 
between the particles. For the marginal case of mass 
independent transfer rates (7 = 0) the interaction van- 
ishes and the exclusion process reduces to coalescing ran- 
dom walkers. This is a well-known problem which can be 
solved exactly e.g. with the method of interparticle dis- 
tribution functions |15fl . The most relevant results for our 
case are: (i) the average ripple mass grows asymptotically 
as (m)(t) y/t and (ii) the ripple size distribution (the 
probability of finding a ripple of mass m at time t) can 
be written in a scaling form as p(m; t) = m G(m/ (m)), 



where the scaling function G(x) = \ixx 2 e KX 



/4 



A length conserving mapping is obtained by consider- 
ing the model with ripple extinction but without removal 
of empty lattice sites. Particles hop over these sites in 
order not to allow for creation of new ripples. Again the 
empty sites become particles in the exclusion picture but 
now the mass transfer corresponds to a hop of a particle 
over all the particles in the same cluster. This is illus- 
trated in Fig. ||. 

A long one particle hop can be considered as moving 
a cluster of particles as a whole. The loss of a ripple 
turns into the aggregation of clusters. In this way our 
model further maps to a cluster-cluster aggregation pro- 
cess where each cluster moves with a rate that depends 
on the distance to its neighbor (s). As the ripples map 
to holes, the main interest lies in the size distribution of 
distances between clusters and not in the cluster size dis- 
tribution itself, as is usually the case. Onc-dimcnsional 
cluster-cluster aggregation models generally obey univer- 
sal dynamical scaling (see Jl(| and references therein), 
which will be seen to be the case also for the sand ripple 
model. 



III. NOISE-INDUCED COARSENING (7 < 0) 

It is known from the mean-field analysis of S that 
for 7 = — 1 the average ripple size grows as (m)(t) « 



ln(t- 



,(m)(Q)-i 



Intuitively, the slow growth follows from 



the fact that, for 7 < 0, the ripples near extinction are 
those with the highest incoming mass rates. Therefore 



the disappearance of a ripple involves a rare fluctuation; 
within the approach of |3| , the mass of a ripple evolving 
in a background of mean mass (m) (t) performs a random 
walk that is biased away from zero. 

We base our theoretical analysis on this observation. 
In what follows, we will assume that at long times, i.e., 
for large (m), the extinction of a ripple is such a rare 
event, that it does not affect the ripple size distribu- 
tion. Neglecting extinctions, we show that the steady 
state distribution is given by a product form. After vali- 
dating the quasi-static approximation by simulations, we 
use it to show that to leading order in t the average rip- 
ple mass grows as (to) ~ — 7 _1 ln(<). Finally we consider 
the approach to the product measure by studying nearest 
neighbor correlations. 



A. Mass distribution 

Without extinctions the ripple size distribution can be 
solved exactly. This is due to the short range of inter- 
actions between ripples: The mass transfer rate depends 
only on t he mass at the target site. As was noted above in 
Sec. II B, this is similar to the zero-range process, where 
the rate depends only on the site of departure. The most 
important characteristics of a zero range process is that 
its steady state is described by a product measure |§. 
This was shown to generalize to processes where the tran- 
sition rate is a product of functions of the occupation 
numbers at the site of departure and the target site |l7j . 
As our model is a special case of this class of models, the 
results of jl7j apply here as well. For completeness we 
give a brief derivation. 

The product measure property implies that the sta- 
tionary probability distribution, P({m;}), of finding the 
system in configuration {mi,TO2, ■ ■ ■ ,tojv} factorizes as 



P({mi}) = Y[p(mi), 



(3) 



where p{rrii) is the probability of finding mass to, at site 
i. In the steady state there are no correlations between 
the ripple sizes. Starting from the master equation for 
P({rrii}) and using the product form (||), one obtains for 
the asymmetric case (the calculation is not presented here 
since up to index changes it is identical to that presented 
in |) 



p(m i )p(m i+ i)T(m l ) = p{m l + l)p(m 



i+i 



l)r(m l+1 - 1). 

(4) 



The condition (Q) , known as pairwise balance jOI , gen- 
eralizes the detailed balance condition familiar from equi- 
librium statistical mechanics. It has a simple interpreta- 
tion. The left hand size of Eq. ([!]) represents the mass 
transfer to the site i which has to balanced by a transfer 
out of this site (the right hand side) in order to be in the 
steady state. The first two terms give the probability to 
find a mass rxii at site i with a right neighbor with mass 
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rrii+i and the last term describes the rate at which the 
site i gains mass from its neighbor. We emphasize that, 
provided a solution to (|4|) can be found, this proves that 
the product measure (|3|) is an exact stationary solution of 
the master equation; on the basis of general arguments, 
this solution is then also expected to be unique. 

Proceeding similarly for the symmetric dynam- 
ics gives (transitions {. . . , mj-i, rtii, J7ii+i, . . . } — > 
{... , m,_i, mj+l,m* + i — 1, . . . } and {. . . , m,_i -1, m,-+ 
l,m i+ i, ...}->•{... , to^i, m n TO i+1 , . . . }) 



p(m 



p(m i _i)Tj(m l )p(mi + i)r(m l ) = 
l)p(mj + l)p(m i+ i)r(m J _i - 1) 



(5) 



Since p(mi+i) cancels out we end up with Eq. (|J). There- 
fore the steady state distribution is independent of the 
asymmetry of the dynamics. 
Equation ([|) can be recast as 

_ x p{mi)T(mi) p(m i+ i - l)T(m i+ i - 1) . . 

a := ; rr— = ; ; , (OJ 



p(m t + 1) 



P(m i+ i) 



where a must be a constant. Denoting p(0) = po and 
recursively iterating equation ([j]), we obtain 

m— 1 

p(m) = Po a m J] r(i) = p a m [(m - l)!] 7 , (7) 



where the product for m = 1 is defined to give unity and 
the last form follows from the definition T(m) = m 7 . 

The unknown constants po and a can be determined by 
the normalization X)m=oP( m ) = 1 ano - ^ ne expectation 
value (to) := X)m=o m P( m )- Explicit results for 7 = 
— 1 and —2 can be found in Appendix |A| for 7 = — 1, 
(fjj) is a (shifted) Poisson distribution. In general, the 
distribution for (771) 3> 1 can be written as 

p(m) = C 2 (7)e 7<m) (m)- 7m -( 1 -T)/ 2 [(m- l)!] 7 , (8) 

where the explicit form of (^2(7) is not important for our 
purposes. Using the form g iven in Eq. (p| ), it is easy to 
show that the width a := yj (m 2 ) — (m) 2 of the distribu- 
tion behaves as a ~ ^/ (m) independent of 7. 

The calculated distributions are compared to numerics 
in Figs. ^ and |5|. The average ripple mass is not a con- 
stant as the simulations include also ripple extinction. 
The excellent agreement at long times shows that indeed 
these become so rare, that between subsequent extinc- 
tions the system has time to equilibrate to the steady 
state. Note that all initial distributions converge to the 
universal distribution p(m;t) given by Eq. (j^), where 
the time-dependence enters only through the mean ripple 
mass (m)(t). 



B. Coarsening law 

Next we proceed to calculate the mean ripple size 
(m)(t) using an approach similar to the analysis of the 
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FIG. 4: The ripple size distributions obtained from simula- 
tions for 7 = -1 at t = 4 (x), 256 (©), 16384 (V), and 
2097152_fD) tog ethe r with the analytical result (solid lines) 
[Eqs. (A2) and (A3)]. The initial distribution at t = is 
a random one and simulations are averaged over 2000 runs 
for a system of size M — 50000. The dashed line shows the 
asymptotic solution given by Eq. (H) for t = 2097152. 




FIG. 5: The ripples size distributions obtained from simu- 
lations for 7 = -0.5 at t = 4 (x), 128 (o), 4096 (V), and 
131072 (□) together with the analytical result (solid lines) 
[Eq. (0)]. The initial distribution at t = is a random one 
and simulations are averaged over 500 runs for a system of size 
M — 50000. The dashed line shows the asymptotic solution 
given by Eq. (§) for t = 131072. 



backgammon model |0| . We assume that, at long times, 
the probability for a given ripple to vanish is equal to 
the probability p(0) obtained by extrapolating the steady 
state probability distribution m) to 777 = 0. The number 
TV of ripples then decays according to dN/dt w — p(0)N. 
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FIG. 6: The growth of the average ripple mass a function of 
time for 7 = -0.5 (*), -0.75 (o), -1 (A), -1.5 (O), and 
—2 (□). The least squares fits are shown by solid lines. The 
inset compares the fitted prefactors (O) to the analytic re- 
sult — 7~ x (solid line). The system sizes range from 50000 to 
100000 and averages are taken over at least 50 independent 
runs. 



Since (m)(t) = M/N we obtain 



w p(0)(m)(t) ~ e^ m >(mr (1 ~ 7)/2 , 



dt 

which to leading order in t gives 

(mXtJw-T-'M*)- 



(9) 



(10) 



Simulations with different initial conditions are in accord 
with Eq. © (Fig. §). 



C. Decay of correlations 

The product measure for the ripple size distribution 
implies that there are no correlations between neighbor- 
ing ripples. This is true only asymptotically. To study 
the approach to the product measure distribution we con- 
sider the normalized nearest neighbor time correlation 
function 



9(t) :-- 



H+l/ 



(11) 



As is clear from Fig. Q the early time behavior is sen- 
sitive to the details of the initial distribution. In this 
regime it is possible to have positive correlation between 
neighboring ripples but at long times there will always 
be anticorrelations, i.e., g(t) < 0. The numerically ob- 
served correlations seem to be independent of the initial 
conditions and vanish in a universal manner as 



/(*) 



-1/2 



FIG. 7: The nearest neighbor correlation function g(t) for 
7 = — 1. The initial condition is random ((m)(0) = 2,D; 
(m)(0) = 1.2,1), monodisperse (rre(0) = 5,o) or Poisson dis- 
tributed «m)(0) = 5,0; (m){0) = 10, A). Open (filled) sym- 
bols correspond to asymmetric (symmetric) dynamics. The 
inset shows the decay at late times for the random case. The 
solid and dashed lines are guides to the eye with slopes —1/2 
and —2/3, respectively. 



for both symmetric and asymmetric dynamics. 

At first sight one may be tempted to relate the de- 
cay of correlations to the extinction events, which per- 
turb the product measure. However, as was shown in 
Sec. [II B , the probability of extinction events decays as 
p(0) ~ e 7 ( m ^*) ~ which is much faster than the nu- 
merically observed decay law (|l2|) . This implies that the 
power law (fl^ ) is associated with the dynamics between 
extinction events, which can be described using standard 
hydrodynamic fluctuation theory for a one-dimensional 
system with a single conserved density. 

Let 4>{x,t) denote the coarse grained mass fluctuations 
in the (quasi-) steady state of mean mass (m) . The long 
wavelength behavior of </> is governed by a Langevin equa- 
tion of the generic form 



dt 



d 2 4> ^dcj) drj 
dx 2 dx dx ' 



(13) 



(12) 



where r}(x, t) is Gaussian white noise with covariance 
{rj(x,t)ri{x',t')) = D5{x-x')S(t_t'). For// ^ 0, Eq. @) 
is the noisy Burgers equation p0[ , which has been widely 
studied in the context of driven diffusive systems |l9f and 
interface growth [^l], ^2], p3| . 

The coefficients v, (i, and D appearing in the long 
wavelength description can be related to the microscopic 
dynamics of the sand ripple model as follows. The non- 
linear term on the right hand side of Eq. (|l^) is gener- 
ated by the asymmetry, and its coefficient is given by 
(X = j"({m)), where j is the steady state mass current. 
Since in our case j = T, we conclude that // ~ (m) 7 ~ 2 . 
In the symmetric case /i = and the diffusion coefficient 



7 



v is proportional to V ~ (m) 7 " 1 (this can be seen by 
expanding Eq. (|l|) around the homogeneous state). Fi- 
nally, owing to a fluctuation-dissipation theorem pfl] , the 
equal time correlations of (13) are Gaussian with covari- 
ance (0(x)0(x')) ~ {D /v)5{x — x') i ndependent of fi. As 
we have shown above in Sec. Ill A , in the ripple model 
the variance of the mass fluctuations is always of order 
(to), hence D/v ~ (to). 

We want to use equation ([l3]) to describe the approach 
to the steady state, starting from some initial condition 
(e.g., the monodisperse state = 0) specified at t = 0. 
The analysis of Eq. (|1|) shows that at long times, and 
for x ^ x', the pair correlation function takes the scaling 
form d| || 

(0(x,f)0(x',t)) = ^^g(^)- 1 |x-x'|). (14) 

Here <J is a scaling function, and £(t) denotes the dynamic 
correlation length. The prefactor of the scaling function 
on the right hand side of Eq. ( O ) is fixed by the require- 
ments that (i) the steady state density fluctuations are 
proportional to D/v, and (ii) the integral over the pair 
correlation function is constant due to mass conservation. 
The correlation length grows diffusively as £(t) ~ (vt) 1 ' 2 
for jit = and superdiffusively as £(i) ~ [(D /v) 1 / 2 [it] 2 / 3 
for /i ^ 0. 

Keeping |x — x'| fixed and taking t — > oo, we see that 
the pair correlations ( |l4| ) decay as (D/v)G(0)^(t)~ 1 . Ex- 
pressing v and D /v in terms of the mean ripple mass, we 
conclude that in the symmetric case (/i = 0) the normal- 
ized correlation function ( |l2] ) should decay as 



)-(l+7)/2 



(lni) 



-(l+T)/2 



t i/a 



(15) 



For 7 = — 1 the logarithmic factor disappears and (ITS) 
becomes a pure power law with exponent —1/2, in accor- 
dance with the simulation results shown in Fig. |?]. More- 
over, the explicit calculation for the diffusive case shows 
that the scaling function Q in Eq. (|lj) is negative, hence 
g{t) < as is observed numerically. 

In the asymmetric case the fluctuation theory pre- 
dicts an asymptotic decay as g(t) ~ l/£(i) ~ t~ 2 / 3 , 
with logarithmic corrections due to the growth of (m) (t) . 
Then why do we find g(t) ~ i^ 1 / 2 also for asymmet- 
ric dynamics? To answer this question, we recall that 
the asymptotic, superdiffusive behavior predicted by the 
noisy Burgers equation sets in only beyond a crossover 
time scale t x , which increases rapidly with decreasing 
strength /1 of the nonlinear term p3| . The crossover 
time is of the order of t x ~ v 5 /{D 1 ^ 4 ). Inserting the 
estimates for /x, v, and D/v derived above, we see that 
for the ripple model 



\3-7 



(In*) 



3-7 



(16) 



For the case 7 = — 1 considered in Fig. |?], this implies that 
superdiffusive behavior can be expected only for times 



such that t/(\nt) 4 ^> 1. The left hand side of this in- 
equality becomes equal to unity for t ss 5500 and reaches 
the value 10 only for t ss 235000. Thus the asymptotic 
regime has not been reached in our simulations. The 
slight deviation of the simulation data from the i -1 / 2 - 
behavior seen after t = 10 5 may indicate the beginning 
of the crossover. 



IV. UNSTABLE COARSENING (7 > 0) 

For 7 > the homogeneous state is linearly unsta- 
ble because the largest ripples are those with the high- 
est growth rate. Ripple extinction is then no longer a 
rare event , and the product measure solution derived in 
Sec. Ill A| becomes invalid. On the other hand, it is plau- 
sible (and will be confirmed by simulations, see below) 
that the linear instability supersedes the noise in the time 
evolution, so that the deterministic equations (|l|) and (||) 
and the stochastic ripple model show the same behavior. 

In what follows, we first develop a mean- field theory 
for the deterministic model in the simplest case of a linear 
robber function (7 = 1). Simulations show that the mean 
field theory is not quantitatively correct, presumably due 
to the neglect of spatial fluctuations. In the nonlinear 
regime < 7 < 1 we use scaling analysis to derive the 
coarsening law. 



A. Mean-field analysis for 7 = 1 

We start our analysis from the deterministic equa- 
tions (Q) and (||). For 7=1 these become linear but 
the system is still non-trivial due to the ripple extinction. 
As the system is deterministic, the only randomness lies 
in the initial condition. We denote the initial ripple size 
distribution by Po(<V>) and its average by Ao. 

The mean field approximation consists of replacing the 
ripples surrounding an arbitrary ripple of size A by ripples 
of the average size (A) , such that the evolution equation 
becomes 



dA 
"dT 



r(A)-r((A» = A-(A). 



(17) 



On this level there is no difference between symmetric 
and asymmetric mass transfer. The solution of Eq. ( |T^ ) 
reads A(Ao,t) = e*[Ao — F(t)], where the function 



F(t) := / dr e-(A)(r) 



(18) 



has to be calculated self-consistently. Note that at this 
point we do not explicitly restrict A(i) to be nonnegative 
(this constraint will enter later). Once F(t) is known, 
the ripple size distribution at time t can be obtained by 
inverting the solution for A(Ao) and inserting this into 
the initial distribution, with the result 



p(X-t) ^e^Poie^X + Fit)). 



(19) 
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Thus in the mean-field approximation the ripple size dis- 
tribution preserves its initial shape but gets scaled and 
shifted. 

It is possible to derive a differential equation for the 
unknown function F(t). The fraction p(t) of surviving 
ripples is equal to the probability that A(t) > 0, 



Pit) 



dA p(A;t) - / dxP (x) =: Pg(F(t)), 

ii JF(t) 



(20) 



where the last equation defines the cumulative distribu- 
tion Pq. The average ripple size is given by (A)(t) = 
Xo/p(t). Inserting this into the definition of F(t) and 
differentiating once gives 



dt 



(21) 



Hence the problem reduces to solving the differential 
equation for a given initial distribution -Po(A). 

For example, for an exponential initial distribution 
P (A ) = XV 1 e~ Ao/Xo , we find F(t) = X t and (X)(t) = 
Aoe*, whereas for a flat distribution 



Po(A ) = 




Ao < 2Ao 
otherwise 



(22) 



the solution is given by F(t) — 2Ao(l — e~*/ 2 ) and 
(X)(t) = Aoe*/ 2 . As the rate of exponential growth is dif- 
ferent in these two cases, we conclude that the coarsening 
behavior of the mean field model ( |l7|) is nonuniversal. 

In general, the exponential growth rate of the mean 
ripple size is governed by the extremal statistics of the ini- 
tial distribution Pq . If the initial ripple sizes are bounded 
by a maximal size A ma x, and Po(Ao) ~ (A max - A ) a for 
Ao — * A max , then the analysis of Eq. ( ^l| ) shows that 
t^ 1 ln(A)(i) (a+ l)/(a + 2), while for fat initial distri- 



An 



we find 



butions with a power law tail, Po(Ao) 
t- 1 ln(A)(t) -> b/{b-l). 

To compare the predictions of the mean-field theory 
to simulations we prefer to show the complement of the 
cumulative distribution 



i(X;t) 



dx p{x;t) =: f 



A 



(A)(t) 



(23) 



where the last equation defines the scaling function f(x). 
A similar definition applies to p(m; t) with the integral 
replaced by a sum. In the case of an exponential distribu- 
tion p(X; t) also the function f(x) is exponential whereas 
for a flat p(X; t) it is linear. 

We solved the deterministic equations ([j]) and (|l7|) us- 
ing the fourth-order Runge-Kutta method |24|j . As a 
check of the algorithm, we reproduced the solution ( |l9| ) 
of the mean-field equations. For the full noiseless sys- 
tem (Q) the exponential initial distribution remains un- 
changed [Fig. ^ a); dashed lines] but also a flat initial 
distribution presumably approaches the exponential one 



m/<m> or X/<k> 
5 10 



L5 




5 10 
m/<m> or )J<X> 



15 



FIG. 8: The complements of the cumulative ripple size distri- 
bution for 7 = 1. a) The distributions for random (exponen- 
tial) initial distribution with noisy (deterministic) dynamics 
are denoted by solid lines (dashed) lines, b) The distributions 
for monodisperse (flat) initial distribution with noisy (deter- 
ministic) dynamics are denoted by solid lines (dashed) lines. 
The curves are shown at times t = 1, . . . ,9 and the thick solid 
lines in both figures represent the function e^ m ^ m K 



[Fig. U b); dashed lines]. This is in conflict with the 
mean-field prediction. In both cases the average ripple 
size grows as (m)(t) ~ e . 



Similarly, in the discrete, noisy ripple model the ran- 
dom initial distribution quickly converges towards an 
exponential scaling function [Fig. || a); solid lines]. 
The monodisperse initial condition spreads out and ap- 
proaches the same form [Fig. || b); solid lines]. Again, 



the mean ripple size grows as (m)(t) 
distributions. 



for both initial 



Since the deterministic model behaves in a similar 
manner as the noisy one, we conclude that, in contrast to 
the case 7 < 0, the noise is irrelevant. The discrepancy 
between the mean-field theory and the full deterministic 
system suggests that the spatial fluctuations are impor- 
tant, as is often the case for low dimensional systems. 
In particular, the numerical results indicate that, in con- 
trast to the mean field prediction, the behavior of the full 
system is universal with respect to the initial ripple size 
distribution. 
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B. 



Coarsening law for < 7 < 1 



As the mean-field equation is not readily solvable for 
7 / 1 and probably would not describe the problem 
correctly anyway, here we present a simple scaling ar- 
gument for the growth of the mean ripple size in the 
regime < 7 < 1. We start from the observation that in 
the linearly unstable case (r'(A) > 0) predominantly ev- 
ery second ripple grows and every second one shrinks j| . 
Therefore we may consider a simplified system consisting 
of two ripples of initial sizes A" > A". We calculate the 
time t* at which the average size has doubled. It is given 
by the conditions A x (t*) = A? + A° and A 2 (t*) = 0. Since 
the mass is conserved we have A :— Xi (t) + X 2 (t) = const.. 

Applying Eq. M gives 




AT, 



(24) 



where the dot denotes derivative with respect to time. 
The solution is implicitly given by \2{t) — A — X\(t) and 



Al « dx 



afi - (A - xp ' 



(25) 



which together with the definition of t* implies the ho- 
mogeneity relation 



t*{a\\,a\) =a 1 -^*(A 1 ,A). 



(26) 



Assuming that the evolving ripple size distribution is gov- 
erned by a single size scale, it follows that the doubling 
time depends on the mean ripple size as t* ~ (A) 1-7 . 
The inverse of the doubling time is the growth rate of 
(A). Hence we may write 



(27) 



which yields (A)(t) ~ t* with z = 1/(1-7). This is 
confirmed by simulations, which give z — 1.32 ±0.02 and 
1.98 ± 0.03 for 7 = 0.25 and 0.50, respectively. We also 
numerically checked the universal scaling behavior of the 
ripple size distribution for < 7 < 1, but in this region 
the scaling function is more complicated than a simple 
exponential. 



V. CONCLUSIONS 

In this paper we have studied a one-dimensional model 
for sand ripple evolution, where mass is transferred be- 
tween neighboring sites with algebraic rates and sites 
containing no mass are removed from the system (rip- 
ple extinction). As the rates depend only on the site to 
which the mass is transferred, the system is similar to 
a zero range process. Thus the steady state in the ab- 
sence of ripple extinction is characterized by a product 



measure. Asymptotically this continues to hold for al- 
gebraically decaying mass transfer rates (7 < 0) since 
the extinctions are exponentially rare at late times. As 
a consequence the average ripple size grows to leading 
order logarithmically slowly with the prefactor — 1/7. 

For 7 < the approach to the steady state product 
measure is algebraic. The correlations between masses of 
neighboring ripples decay universally, i.e., independent 
of the initial distribution, as t -1 / 2 and i -2 / 3 for sym- 
metric and asymmetric mass transfer, respectively. In 
the asymmetric case the asymptotic regime is preceded 
by a long crossover where i _1 / 2 -decay is observed. 

For algebraically growing robber functions (7 > 0) the 
coarsening is driven by the linear instability of the ho- 
mogeneous state. Ripple extinctions become frequent, 
and the product measure is no longer relevant. The 
average ripple size grows algebraically as for 
< 7 < 1. The behavior at 7 = is discontinuous since 
(m)(t) ~ i 1 / 2 for 7 = 0, which follows from the mapping 
to coalescing random walkers. For 7=1 (m) (i) ~ e and 
the scaling function of the ripple size distribution appears 
to be a simple exponential. The dynamical noise, which 
is necessary to have coarsening for 7 < 0, is irrelevant 
for 7 > 0. The mean- field theory developed for 7=1 
reproduces the exponential growth of the mean ripple 
size, but it is insufficient to describe the universality of 
the growth law and the ripple size distribution which is 
observed numerically. 

It is interesting to compare the results to the behav- 
ior in one-dimensional cluster-cluster aggregation. Re- 
call that the model treated here can be mapped to 
cluster-cluster aggregation with hopping rates of clus- 
ters depending algebraically on the distance between 
them 



(Sec, 
m 7 



1IB). When the hopping rates depend as 



T(m) ~ m 7 on the masses of clusters, the growth of the 
average cluster size is algebraic with (m)(t) ~ t 1 /( 2 -i) 
for all 7 < 2 [jl6). Thus the behavior for nonnegative 
values of 7 is rather similar in the two models, but for 
7 < one finds a drastic difference due to the repulsive 
interaction between clusters in the ripple model. 

We conclude by adducing some open problems for fu- 
ture studies. One of the most interesting issues is to 
understand the coarsening and the final ripple size se- 
lection in the case of a nonmonotonic robber function. 
This has direct applications in the coarsening of vor- 
tex ripples, where the robber function has recently been 
measured j2j. Initially these systems are in the unsta- 
ble regime, where the transfer function is monotonically 
increasing. As we have seen in the present article, even 
this is a harder problem than the case where the domi- 
nant contribution to coarsening comes from the dynami- 
cal fluctuations. For nonmonotonic robber functions one 
needs an understanding of both coarsening mechanisms. 
Therefore the starting point into this direction would be 
to better understand the 7 > case. 
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unity, and we set (— 1)! — 1. The normalization condition 

Em=oP( m ) = 1 S iveS 



pirn) 



(1 + ae a )- 1 a m /(m- 1)! for 7 = -1 
I ~ 1 (2V^)a m ~ 1 /[(m - l)!] 2 for 7 = -2, 



(A2) 



where I n (x) is the modified Bessel function of the first 
kind. The parameter a is related to the expectation value 

(™>) = Em=o rnp(m) by 



\a(a+l)e a /(l + ae°) for 7 = -1 (A ) 
y ' I 1 + y/Sh{2y/a)/Io{2y/a) for 7 = -2. ' 



APPENDIX A: CALCULATION OF MASS 
DISTRIBUTIONS FOR 7=1 AND 2 

Here we calculate the explicit form of the mass distri- 
bution 



p(m) = Po a m 11 r(i)= Po a m [(m-l)\y, (Al) 

i=l 

for T(m) = m? in the case s 7 = —1 and —2. For m = 
1 the product in Eq. ( |AlD for m = 1 is defined to be 



Using the expansions I n (x) = e x / \/2irx+0(l/x) for a — > 
00, these formulae simplify to (m) « a and (m) ~ ^/a, 
for 7 = — 1 and —2, respectively. Hence for (m) — > 00 
the distributions become 



pirn) 



-< m >(TO) m -V(m - 1)! 



for 7 = — 1 



2 v /^e- 2 < m >(™) 2m - 3 /7[(m- l)!] 2 for 7 



(A4) 



which are of the general form indicated in Eq. m) . 
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